library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ readr 2.1.4
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.4.4 ✔ tibble 3.2.1
## ✔ lubridate 1.9.3 ✔ tidyr 1.3.0
## ✔ purrr 1.0.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(ggpubr)
library(rstatix)
##
## Attaching package: 'rstatix'
##
## The following object is masked from 'package:stats':
##
## filter
library(nparLD)
## Loading required package: MASS
##
## Attaching package: 'MASS'
##
## The following object is masked from 'package:rstatix':
##
## select
##
## The following object is masked from 'package:dplyr':
##
## select
library(PMCMRplus)
library(WRS2)
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib
import seaborn as sns
import numpy as np
from scipy.stats import pearsonr, spearmanr
import os
import gseapy as gp
from gseapy.parser import gsea_gmt_parser
from scipy.stats import wilcoxon, ttest_ind, ttest_rel
from statsmodels.stats.multitest import multipletests
from glob import glob
from statannotations.Annotator import Annotator
plt.rcParams.update({'font.size': 14})
annot = pd.read_csv('../data/batches12_chemo_regimen_annotation.tsv', sep='\t', index_col=0)
annot['Primary Study ID'] = annot.loc[:, 'Primary Study ID'].astype(str)
pt_reg = annot.set_index('Primary Study ID').loc[:, 'chemo_regimen'].to_dict()
platinum_type_dict = annot.set_index('Primary Study ID').loc[:, 'platinum_type'].to_dict()
def dyn_comp(u, ax, lin_dens3, log_dens3, title, xlabel, test='t-test_paired', test_name = 't-test paired', size=3, format_simple = True):
sns.boxplot(y=u, x='sample_type', data = lin_dens3, ax=ax, order = ['TUR', 'Cystectomy'],
palette = tur_cys_palette, showfliers=False)
sns.stripplot(y=u, x='sample_type', data = lin_dens3, ax=ax, order = ['TUR', 'Cystectomy'],
color='black', size=size)
sns.lineplot(y=u, x='sample_type', data = lin_dens3, ax=ax, units = 'pt_number', estimator=None, linewidth=0.35, color = 'grey')
ax.set_yscale("log")
if format_simple:
ax.yaxis.set_major_formatter(matplotlib.ticker.ScalarFormatter())
ax.set_xlabel('')
ax.set_ylabel(xlabel)
ax.set_title(title)
annotator = Annotator(ax, pairs = [('TUR', 'Cystectomy')], data = log_dens3, x='sample_type', y=u, order=['TUR', 'Cystectomy'])
annotator.configure(test=test, text_format = 'simple', comparisons_correction='BH')
_, test_results = annotator.apply_test()._get_output()
pv = test_results[0].data.pvalue
annotator.set_custom_annotations(['{} p={:.2g}\nN = {}'.format(test_name, pv, len(log_dens3.index))])
annotator.annotate()
print(u, 'Median TUR = {:.3g}'.format(lin_dens3[lin_dens3.sample_type == 'TUR'].loc[:, u].median()),
'\nMedian Cystectomy = {:.3g}'.format(lin_dens3[lin_dens3.sample_type == 'Cystectomy'].loc[:, u].median()))
def dyn_comp_lin(u, ax, lin_dens3, title, xlabel, test='t-test_paired', test_name = 't-test paired', size=3):
sns.boxplot(y=u, x='sample_type', data = lin_dens3, ax=ax, order = ['TUR', 'Cystectomy'],
palette = tur_cys_palette, showfliers=False)
sns.stripplot(y=u, x='sample_type', data = lin_dens3, ax=ax, order = ['TUR', 'Cystectomy'],
color='black', size=size)
sns.lineplot(y=u, x='sample_type', data = lin_dens3, ax=ax, units = 'pt_number', estimator=None, linewidth=0.35, color = 'grey')
ax.set_xlabel('')
ax.set_ylabel(xlabel)
ax.set_title(title)
annotator = Annotator(ax, pairs = [('TUR', 'Cystectomy')], data = lin_dens3, x='sample_type', y=u, order=['TUR', 'Cystectomy'])
annotator.configure(test=test, text_format = 'simple', comparisons_correction='BH')
_, test_results = annotator.apply_test()._get_output()
pv = test_results[0].data.pvalue
annotator.set_custom_annotations(['{} p={:.2g}\nN = {}'.format(test_name, pv, len(lin_dens3.index))])
annotator.annotate()
print(u, 'Median TUR = {:.3g}'.format(lin_dens3[lin_dens3.sample_type == 'TUR'].loc[:, u].median()),
'\nMedian Cystectomy = {:.3g}'.format(lin_dens3[lin_dens3.sample_type == 'Cystectomy'].loc[:, u].median()))
def dyn_comp_delta(u, ax, lin_dens3, title, xlabel, to_compare, order, test = 't-test_paired', size=3, ofs = 1):
sns.boxplot(y=u, x=to_compare, data = lin_dens3, ax=ax, order = order,
color = 'mediumseagreen', showfliers=False)
sns.swarmplot(y=u, x=to_compare, data = lin_dens3, ax=ax, order = order,
color='black', size=size)
sns.lineplot(y=u, x=to_compare, data = lin_dens3, ax=ax, units = 'pt_number', estimator=None, linewidth=0.7)
ax.set_xlabel('')
ax.set_ylabel(xlabel)
ax.set_title(title)
annotator = Annotator(ax, pairs = [order], data = lin_dens3, x=to_compare, y=u, order=order)
annotator.configure(test=test, text_format = 'simple', comparisons_correction='BH', line_offset = ofs, line_offset_to_group = ofs)
_, test_results = annotator.apply_test()._get_output()
pv = test_results[0].data.pvalue
annotator.set_custom_annotations(['{} p={:.2g}'.format(test, pv)])
annotator.annotate(line_offset = ofs, line_offset_to_group = ofs)
norm_log_dens = pd.read_csv('../data/combined_norm_log10_density.tsv', sep='\t', index_col=0)
norm_log_dens_tumor = pd.read_csv('../data/tumor_norm_log10_density.tsv', sep='\t', index_col=0)
norm_log_dens_stroma = pd.read_csv('../data/stroma_norm_log10_density.tsv', sep='\t', index_col=0)
raw_dens = pd.read_csv('../data/all_densities_combined_kde_mm2.tsv', sep='\t', index_col=0)
dens_to_comp = ['B_cells', 'CD4_Tcells', 'Macrophages', 'Tregs', 'CD8_Tcells', 'PanCK+_cells', 'Negative_cells']
def transform_dens(norm_log_dens, addition = 0):
norm_log_dens = norm_log_dens.loc[:, dens_to_comp] + addition
lin_dens2 = (10**norm_log_dens).assign(sample_type = norm_log_dens.index.str.split(' ').str[1]).assign(pt_number = norm_log_dens.index.str.split(' ').str[0])
lin_dens3 = lin_dens2[lin_dens2.pt_number.isin(lin_dens2.pt_number.value_counts()[lin_dens2.pt_number.value_counts() == 2].index)].sort_values(['pt_number', 'sample_type'], ascending=False)
lin_dens3 = lin_dens3.assign(chemo_regimen = lin_dens3.pt_number.map(pt_reg)).assign(platinum_type = lin_dens3.pt_number.map(platinum_type_dict))
log_dens2 = norm_log_dens.assign(sample_type = norm_log_dens.index.str.split(' ').str[1]).assign(pt_number = norm_log_dens.index.str.split(' ').str[0])
log_dens3 = log_dens2[log_dens2.pt_number.isin(log_dens2.pt_number.value_counts()[log_dens2.pt_number.value_counts() == 2].index)].sort_values(['pt_number', 'sample_type'], ascending=False)
log_dens3 = log_dens3.assign(chemo_regimen = log_dens3.pt_number.map(pt_reg)).assign(platinum_type = log_dens3.pt_number.map(platinum_type_dict))
return lin_dens3, log_dens3
lin_dens3, log_dens3 = transform_dens(norm_log_dens, addition = 2)
raw_lin_dens3, raw_log_dens3 = transform_dens(raw_dens)
lin_dens3_tumor, log_dens3_tumor = transform_dens(norm_log_dens_tumor, addition = 2)
lin_dens3_stroma, log_dens3_stroma = transform_dens(norm_log_dens_stroma, addition = 2)
pdl1 = pd.read_csv('../data/Prepost NAC patients overview PD-L1.csv', sep=';', index_col=0)
pdl1_cys = pdl1.loc[~pdl1.index.isna(), ['Material.1', 'CPS (%).1', 'IC (%).1', 'TC (%).1']].dropna()
pdl1_tur = pdl1.loc[~pdl1.index.isna(), ['Material', 'CPS (%)', 'IC (%)', 'TC (%)']].dropna().replace('76,5', '76.5')
pdl1_tur.index = pdl1_tur.index.astype(int).astype(str) + ' TUR'
pdl1_cys.index = pdl1_cys.index.astype(int).astype(str) + ' Cystectomy'
pdl1_cys.columns = ['Material', 'CPS (%)', 'IC (%)', 'TC (%)']
pdl1_tur = pdl1_tur.drop('Material', axis=1).astype(float)+1
pdl1_cys = pdl1_cys.drop('Material', axis=1).astype(float)+1
pdl1_df = pd.concat([pdl1_tur.assign(sample_type = 'TUR'), pdl1_cys.assign(sample_type = 'Cystectomy')])
pdl1_df_log = pd.concat([np.log10(pdl1_tur).assign(sample_type = 'TUR'), np.log10(pdl1_cys).assign(sample_type = 'Cystectomy')])
pdl1_df2 = pdl1_df.assign(pt_number = pdl1_df.index.str.split(' ').str[0])
pdl1_df3 = pdl1_df2[pdl1_df2.pt_number.isin(pdl1_df2.pt_number.value_counts()[pdl1_df2.pt_number.value_counts() == 2].index)].sort_values(['pt_number', 'sample_type'], ascending=False)
pdl1_df3 = pdl1_df3.assign(platinum_type = pdl1_df3.pt_number.map(platinum_type_dict)).assign(chemo_regimen = pdl1_df3.pt_number.map(pt_reg))
pdl1_df2_log = pdl1_df_log.assign(pt_number = pdl1_df_log.index.str.split(' ').str[0])
pdl1_df3_log = pdl1_df2_log[pdl1_df2_log.pt_number.isin(pdl1_df2_log.pt_number.value_counts()[pdl1_df2_log.pt_number.value_counts() == 2].index)].sort_values(['pt_number', 'sample_type'], ascending=False)
pdl1_df3_log = pdl1_df3_log.assign(platinum_type = pdl1_df3_log.pt_number.map(platinum_type_dict)).assign(chemo_regimen = pdl1_df3_log.pt_number.map(pt_reg))
ssgsea = pd.read_csv('../data/ssgsea_scores_Fig2.tsv', sep='\t', index_col=0)
ssgsea = ssgsea.reset_index().assign(sample_type = ssgsea.reset_index().loc[:, 'index'].apply(lambda x: x.split(' ')[1])).assign(pt_number = ssgsea.reset_index().loc[:, 'index'].apply(lambda x: x.split(' ')[0])).set_index('index')
ps_ssgsea0 = ssgsea[ssgsea.pt_number.isin(ssgsea.pt_number.value_counts()[ssgsea.pt_number.value_counts() == 2].index)].sort_values(by = ['pt_number', 'sample_type'], ascending=False)
ps_ssgsea0 = ps_ssgsea0.assign(chemo_regimen = ps_ssgsea0.pt_number.map(pt_reg)).assign(platinum_type = ps_ssgsea0.pt_number.map(platinum_type_dict))
tur_cys_palette = {'TUR': 'lightcoral', 'Cystectomy': 'lightskyblue'}
plt.close()
fig, axs = plt.subplots(3, 3, figsize = (15, 15))
af = axs.flat
dyn_comp('CD8_Tcells', af[0], lin_dens3_tumor, log_dens3_tumor, '', 'CD8 T cell percentage in tumor, %')
## TUR vs. Cystectomy: t-test paired p=0.015
## N = 170
## CD8_Tcells Median TUR = 1.33
## Median Cystectomy = 1.95
dyn_comp('CD8_Tcells', af[1], lin_dens3_stroma, log_dens3_stroma, '', 'CD8 T cell percentage in stroma, %')
## TUR vs. Cystectomy: t-test paired p=0.017
## N = 170
## CD8_Tcells Median TUR = 5.1
## Median Cystectomy = 6.51
dyn_comp('CD8_Tcells', af[2], lin_dens3, log_dens3, '', 'CD8 T cell percentage in tumor + stroma, %')
## TUR vs. Cystectomy: t-test paired p=0.00076
## N = 170
## CD8_Tcells Median TUR = 2.34
## Median Cystectomy = 4.17
_ = af[0].set_ylim((0.01), 400)
_ = af[1].set_ylim((0.01), 400)
_ = af[2].set_ylim((0.01), 400)
for i in range(3):
_ = af[i].set_yticks([0.01,0.1,1, 10, 100], ['0.01','0.1','1', '10', '100'])
for i,u in enumerate(['IC (%)', 'TC (%)', 'CPS (%)']):
dyn_comp(u, af[i+3], pdl1_df3, pdl1_df3_log, '', '{} PDL1 score'.format(u), test='Wilcoxon', test_name='Wilcoxon')
_ = af[i+3].set_ylim((0.5, 400))
_ = af[i+3].set_yticks([1,10,100], ['1','10', '100'])
ps_ssgsea = ps_ssgsea0.copy()
dyn_comp_lin('TGFB_Mariathasan', af[6], ps_ssgsea, '', 'fibroblast TGFb ssGSEA score')
dyn_comp_lin('IFNG_signature', af[7], ps_ssgsea, '', 'INFg ssGSEA score')
dyn_comp_lin('CD8_T_EFFECTOR', af[8], ps_ssgsea, '', 'effector CD8 T cell ssGSEA score')
for i in range(9):
_ = af[i].set_xticks(['TUR', 'Cystectomy'], ['TUR-BT', 'Cystectomy'])
plt.subplots_adjust(hspace=0.2, wspace=0.4)
plt.show()
# Figure S3
plt.close()
fig, axs = plt.subplots(2, 6, figsize = (22.7, 11))
af = axs.flat
cell_names = {'CD8_Tcells': 'CD8 T cell', 'Macrophages': 'Macrophage', 'CD4_Tcells': 'T helper', 'Tregs': 'Treg', 'B_cells': 'B cell', 'PanCK+_cells': 'Cancer cell', 'Negative_cells': 'Negative cell'}
for i, u in enumerate(['CD8_Tcells', 'Macrophages', 'CD4_Tcells', 'Tregs', 'B_cells', 'PanCK+_cells']):
dyn_comp(u, af[i], lin_dens3, log_dens3, '', '{} percentage, %'.format(cell_names[u]), test='Wilcoxon', test_name='Wilcoxon')
_ = af[i].set_xticks(['TUR', 'Cystectomy'], ['TUR-BT', 'Cystectomy'])
for i, u in enumerate(['CD8_Tcells', 'Macrophages', 'CD4_Tcells', 'Tregs', 'B_cells', 'PanCK+_cells']):
dyn_comp(u, af[i+6], raw_lin_dens3, raw_log_dens3, '', '{} density, cells/mm2'.format(cell_names[u]), test='Wilcoxon', test_name='Wilcoxon', format_simple = False)
_ = af[i+6].set_xticks(['TUR', 'Cystectomy'], ['TUR-BT', 'Cystectomy'])
plt.subplots_adjust(hspace=0.2, wspace=0.47)
plt.show()
cys_den = log_dens3[log_dens3.sample_type == 'Cystectomy'].set_index('pt_number')
tur_den = log_dens3[log_dens3.sample_type == 'TUR'].set_index('pt_number')
common_in = cys_den.index.intersection(tur_den.index)
delta_den = cys_den.loc[common_in].drop(['sample_type', 'chemo_regimen', 'platinum_type'], axis=1) - tur_den.loc[common_in].drop(['sample_type', 'chemo_regimen', 'platinum_type'], axis=1)
delta_den = delta_den.assign(pt_number = delta_den.index)
delta_den = delta_den.assign(chemo_regimen = delta_den.pt_number.map(pt_reg)).assign(platinum_type = delta_den.pt_number.map(platinum_type_dict))
cys_ps= ps_ssgsea0[ps_ssgsea0.sample_type == 'Cystectomy'].set_index('pt_number')
tur_ps= ps_ssgsea0[ps_ssgsea0.sample_type == 'TUR'].set_index('pt_number')
common_in = cys_ps.index.intersection(tur_ps.index)
delta_ps = cys_ps.loc[common_in].drop(['sample_type', 'chemo_regimen', 'platinum_type'], axis=1) - tur_ps.loc[common_in].drop(['sample_type', 'chemo_regimen', 'platinum_type'], axis=1)
delta_ps = delta_ps.assign(pt_number = delta_ps.index)
delta_ps = delta_ps.assign(chemo_regimen = delta_ps.pt_number.map(pt_reg)).assign(platinum_type = delta_ps.pt_number.map(platinum_type_dict))
all_med = pd.read_csv('../data/all_median_distances_wo_offset.tsv', sep='\t', index_col=0)
all_med_log = np.log10(all_med.drop(['sample_type', 'pt_number'], axis=1)).assign(sample_type = all_med.index.str.split(' ').str[1]).assign(pt_number = all_med.index.str.split(' ').str[0])
all_med_lin = all_med.drop(['sample_type'], axis=1).assign(sample_type = all_med.index.str.split(' ').str[1]).assign(pt_number = all_med.index.str.split(' ').str[0])
all_med_log = all_med_log[all_med_log.pt_number.isin(all_med_log.pt_number.value_counts()[all_med_log.pt_number.value_counts() == 2].index)].sort_values(['pt_number', 'sample_type'], ascending=False).assign(chemo_regimen = all_med_log.pt_number.map(pt_reg)).assign(platinum_type =all_med_log.pt_number.map(platinum_type_dict))
all_med_lin = all_med_lin[all_med_lin.pt_number.isin(all_med_lin.pt_number.value_counts()[all_med_lin.pt_number.value_counts() == 2].index)].sort_values(['pt_number', 'sample_type'], ascending=False).assign(chemo_regimen = all_med_lin.pt_number.map(pt_reg)).assign(platinum_type = all_med_lin.pt_number.map(platinum_type_dict))
plt.close()
fig, axs = plt.subplots(1, 2, figsize = (8, 5))
af = axs.flat
dyn_comp('CD8_Tcell_to_tumor_median', af[0], all_med_lin, all_med_log, '', 'CD8 T cell to 1-NN cancer cell\nmedian distance, um', test='Wilcoxon', test_name='Wilcoxon')
dyn_comp('Macrophage_to_tumor_median', af[1], all_med_lin, all_med_log, '', 'Macrophage to 1-NN cancer cell\nmedian distance, um', test='Wilcoxon', test_name='Wilcoxon')
af[0].set_yticks([10, 50, 100, 500], ['10','50','100','500'])
af[1].set_yticks([10, 50, 100], ['10','50','100'])
#af[0].set_ylim((0.0001), 2)
#af[1].set_ylim((0.0001), 2)
for i in range(2):
af[i].set_xticks(['TUR', 'Cystectomy'], ['TUR-BT', 'Cystectomy'])
#_ = fig.set_tight_layout(True)
plt.subplots_adjust(hspace=0.2, wspace=0.5)
plt.show()
tur_dist2 = (all_med_log[all_med_log.sample_type == 'TUR']).set_index('pt_number')
cys_dist2 = (all_med_log[all_med_log.sample_type == 'Cystectomy']).set_index('pt_number')
common_index0 = tur_dist2.index.intersection(cys_dist2.index)
delta_dist0 = cys_dist2.loc[common_index0, ['CD8_Tcell_to_tumor_median', 'Macrophage_to_tumor_median']] - tur_dist2.loc[common_index0, ['CD8_Tcell_to_tumor_median', 'Macrophage_to_tumor_median']]
delta_dist0 = delta_dist0.assign(pt_number = delta_dist0.index.str.split(' ').str[0])
delta_dist0 = delta_dist0.assign(chemo_regimen = delta_dist0.pt_number.map(pt_reg)).assign(platinum_type = delta_dist0.pt_number.map(platinum_type_dict)).dropna(subset = ['chemo_regimen'])
pdl1_df3_log = pdl1_df3_log.rename(columns = {'IC (%)': 'IC', 'CPS (%)': 'CPS', 'TC (%)': 'TC'}).dropna()
pdl1_df3 = pdl1_df3.rename(columns = {'IC (%)': 'IC', 'CPS (%)': 'CPS', 'TC (%)': 'TC'}).dropna()
cys_pdl = pdl1_df3_log[pdl1_df3_log.sample_type == 'Cystectomy'].set_index('pt_number')
tur_pdl = pdl1_df3_log[pdl1_df3_log.sample_type == 'TUR'].set_index('pt_number')
common_in = cys_pdl.index.intersection(tur_pdl.index)
delta_pdl = cys_pdl.loc[common_in].drop(['sample_type', 'chemo_regimen', 'platinum_type'], axis=1) - tur_pdl.loc[common_in].drop(['sample_type', 'chemo_regimen', 'platinum_type'], axis=1)
delta_pdl = delta_pdl.assign(pt_number = delta_pdl.index)
delta_pdl = delta_pdl.assign(chemo_regimen = delta_pdl.pt_number.map(pt_reg)).assign(platinum_type = delta_pdl.pt_number.map(platinum_type_dict))
def dyn_comp_anova(u, ax, lin_dens3, log_dens3, title, xlabel, loc_leg, test='t-test_paired', test_name='t-test paired'):
sns.boxplot(y=u, x = 'chemo_regimen', hue='sample_type', data = lin_dens3, ax=ax, hue_order = ['TUR', 'Cystectomy'], order = ['MVAC', 'gemcitabine_platinum'],
palette = tur_cys_palette, showfliers=False)
sns.stripplot(y=u, x = 'chemo_regimen', hue='sample_type', data = lin_dens3, ax=ax, hue_order = ['TUR', 'Cystectomy'], order = ['MVAC', 'gemcitabine_platinum'],
color='black', dodge=True, size=3)
#sns.lineplot(y=u, x='chemo_regimen', hue='sample_type', data = lin_dens3, ax=ax, units = 'pt_number', estimator=None, linewidth=0.5, color = 'grey')
ax.set_yscale("log")
ax.set_xlabel('')
ax.set_ylabel(xlabel)
ax.set_title(title)
annotator = Annotator(ax, pairs = [(('MVAC', 'TUR'), ('MVAC', 'Cystectomy')), (('gemcitabine_platinum', 'TUR'), ('gemcitabine_platinum', 'Cystectomy'))], data = log_dens3, x='chemo_regimen', hue = 'sample_type', y=u, hue_order = ['TUR', 'Cystectomy'], order = ['MVAC', 'gemcitabine_platinum'])
annotator.configure(test=test, text_format = 'simple', comparisons_correction='BH')
_, test_results = annotator.apply_test()._get_output()
#pv = test_results[0].data.pvalue
annotator.set_custom_annotations(['{} p={:.2g}'.format(test_name, test_results[0].data.pvalue), '{} p={:.2g}'.format(test_name, test_results[1].data.pvalue)])
annotator.annotate()
handles, labels = ax.get_legend_handles_labels()
ax.legend(handles[0:2], ['TUR-BT', 'Cystectomy'], loc=loc_leg)
def dyn_comp_lin_anova(u, ax, lin_dens3, title, xlabel, loc_leg, test='t-test_paired', test_name='t-test paired'):
sns.boxplot(y=u, x = 'chemo_regimen', hue='sample_type', data = lin_dens3, ax=ax, hue_order = ['TUR', 'Cystectomy'], order = ['MVAC', 'gemcitabine_platinum'],
palette = tur_cys_palette, showfliers=False)
sns.stripplot(y=u, x = 'chemo_regimen', hue='sample_type', data = lin_dens3, ax=ax, hue_order = ['TUR', 'Cystectomy'], order = ['MVAC', 'gemcitabine_platinum'],
color='black', dodge=True, size=3)
#sns.lineplot(y=u, hue='chemo_regimen', x='sample_type', data = lin_dens3, ax=ax, units = 'pt_number', estimator=None, linewidth=0.5, color = 'grey')
ax.set_xlabel('')
ax.set_ylabel(xlabel)
ax.set_title(title)
annotator = Annotator(ax, pairs = [(('MVAC', 'TUR'), ('MVAC', 'Cystectomy')), (('gemcitabine_platinum', 'TUR'), ('gemcitabine_platinum', 'Cystectomy'))], data = lin_dens3, x='chemo_regimen', hue = 'sample_type', y=u, hue_order = ['TUR', 'Cystectomy'], order = ['MVAC', 'gemcitabine_platinum'])
annotator.configure(test=test, text_format = 'simple', comparisons_correction='BH')
_, test_results = annotator.apply_test()._get_output()
#pv = test_results[0].data.pvalue
annotator.set_custom_annotations(['{} p={:.2g}'.format(test_name, test_results[0].data.pvalue), '{} paired p={:.2g}'.format(test_name, test_results[1].data.pvalue)])
annotator.annotate()
handles, labels = ax.get_legend_handles_labels()
ax.legend(handles[0:2], ['TUR-BT', 'Cystectomy'], loc=loc_leg)
plt.rcParams.update({'font.size': 15})
testt = 'Wilcoxon'
testt_name = 'Wilcoxon'
plt.close()
fig, axs = plt.subplots(2, 2, figsize = (14.5, 12.5), width_ratios = [10, 5])
af = axs.flat
u = 'Macrophages'
dyn_comp_anova(u, af[0], lin_dens3, log_dens3, '', 'Macrophage cell percentage', loc_leg=(0.3, 0.03), test='Wilcoxon', test_name = 'Wilcoxon')
dyn_comp_delta(u, af[1], delta_den, '', 'Δ log10(macrophage percentage)\nupon treatment', to_compare = 'chemo_regimen', order = ('MVAC', 'gemcitabine_platinum'), test = 'Mann-Whitney', ofs = -0.08, size=5)
af[0].set_xticklabels(['MVAC', 'gemcitabine +\nplatinum'])
af[1].set_xticklabels(['MVAC', 'gemcitabine +\nplatinum'])
af[1].set_xlim((-0.7,1.7))
af[1].set_ylim((-1.4,2.8))
af[0].set_yticks([0.1,1, 10, 100], ['0.1','1', '10', '100'])
u = 'TGFB_Mariathasan'
dyn_comp_lin_anova(u, af[2], ps_ssgsea0, '', 'Fibroblast TGFb ssGSEA score', loc_leg='lower right', test='Wilcoxon', test_name = 'Wilcoxon')
dyn_comp_delta(u, af[3], delta_ps, '', 'Δ fibroblast TGFb upon treatment', to_compare = 'chemo_regimen', order = ('MVAC', 'gemcitabine_platinum'), test = 'Mann-Whitney', ofs=0.1, size=5)
af[2].set_xticklabels(['MVAC', 'gemcitabine +\nplatinum'])
af[3].set_xticklabels(['MVAC', 'gemcitabine +\nplatinum'])
af[3].set_xlim((-0.7,1.7))
#_ = fig.set_tight_layout(True)
plt.subplots_adjust(hspace=0.3, wspace=0.25)
plt.show()
lin_dens3_cr = lin_dens3[lin_dens3.chemo_regimen.isin(['MVAC', 'gemcitabine_platinum'])]
#lin_dens3_ptt = lin_dens3[lin_dens3.platinum_type.isin(['Cisplatin', 'Carboplatin'])]
log_dens3_cr = log_dens3[log_dens3.chemo_regimen.isin(['MVAC', 'gemcitabine_platinum'])]
log_dens3_ptt = log_dens3[log_dens3.platinum_type.isin(['Cisplatin', 'Carboplatin'])]
py$log_dens3_cr %>%
group_by(sample_type, chemo_regimen) %>%
shapiro_test(CD8_Tcells)
## # A tibble: 4 × 5
## sample_type chemo_regimen variable statistic p
## <chr> <chr> <chr> <dbl> <dbl>
## 1 Cystectomy MVAC CD8_Tcells 0.985 0.936
## 2 Cystectomy gemcitabine_platinum CD8_Tcells 0.984 0.745
## 3 TUR MVAC CD8_Tcells 0.969 0.499
## 4 TUR gemcitabine_platinum CD8_Tcells 0.984 0.759
py$log_dens3_ptt %>%
group_by(sample_type, platinum_type) %>%
shapiro_test(CD8_Tcells)
## # A tibble: 4 × 5
## sample_type platinum_type variable statistic p
## <chr> <chr> <chr> <dbl> <dbl>
## 1 Cystectomy Carboplatin CD8_Tcells 0.910 0.351
## 2 Cystectomy Cisplatin CD8_Tcells 0.991 0.909
## 3 TUR Carboplatin CD8_Tcells 0.923 0.454
## 4 TUR Cisplatin CD8_Tcells 0.981 0.370
py$log_dens3_cr %>%
group_by(sample_type) %>%
levene_test(CD8_Tcells ~ chemo_regimen)
## Warning: There were 2 warnings in `mutate()`.
## The first warning was:
## ℹ In argument: `data = map(.data$data, .f, ...)`.
## Caused by warning in `leveneTest.default()`:
## ! group coerced to factor.
## ℹ Run `dplyr::last_dplyr_warnings()` to see the 1 remaining warning.
## # A tibble: 2 × 5
## sample_type df1 df2 statistic p
## <chr> <int> <int> <dbl> <dbl>
## 1 Cystectomy 1 78 0.296 0.588
## 2 TUR 1 78 1.26 0.264
py$log_dens3_ptt %>%
group_by(sample_type) %>%
levene_test(CD8_Tcells ~ platinum_type)
## Warning: There were 2 warnings in `mutate()`.
## The first warning was:
## ℹ In argument: `data = map(.data$data, .f, ...)`.
## Caused by warning in `leveneTest.default()`:
## ! group coerced to factor.
## ℹ Run `dplyr::last_dplyr_warnings()` to see the 1 remaining warning.
## # A tibble: 2 × 5
## sample_type df1 df2 statistic p
## <chr> <int> <int> <dbl> <dbl>
## 1 Cystectomy 1 74 0.753 0.388
## 2 TUR 1 74 2.88 0.0939
box_m(py$log_dens3_cr[, "CD8_Tcells", drop = FALSE], py$log_dens3_cr$chemo_regimen)
## # A tibble: 1 × 4
## statistic p.value parameter method
## <dbl> <dbl> <dbl> <chr>
## 1 1.10 0.294 1 Box's M-test for Homogeneity of Covariance Matric…
box_m(py$log_dens3_ptt[, "CD8_Tcells", drop = FALSE], py$log_dens3_ptt$chemo_regimen)
## # A tibble: 1 × 4
## statistic p.value parameter method
## <dbl> <dbl> <dbl> <chr>
## 1 1.22 0.270 1 Box's M-test for Homogeneity of Covariance Matric…
res.aov <- anova_test(
data = py$log_dens3_cr, dv = CD8_Tcells, wid = pt_number,
between = chemo_regimen, within = sample_type
)
get_anova_table(res.aov)
## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 chemo_regimen 1 78 0.048 0.827000 0.000431
## 2 sample_type 1 78 13.611 0.000414 * 0.050000
## 3 chemo_regimen:sample_type 1 78 0.857 0.358000 0.003000
bwtrim(CD8_Tcells ~ chemo_regimen*sample_type, id = pt_number, data = py$log_dens3_cr)
## Call:
## bwtrim(formula = CD8_Tcells ~ chemo_regimen * sample_type, id = pt_number,
## data = py$log_dens3_cr)
##
## value df1 df2 p.value
## chemo_regimen 0.0381 1 56.2559 0.8459
## sample_type 8.3784 1 49.8259 0.0056
## chemo_regimen:sample_type 0.2000 1 49.8259 0.6566
res.aov <- anova_test(
data = py$log_dens3_ptt, dv = CD8_Tcells, wid = pt_number,
between = platinum_type, within = sample_type
)
get_anova_table(res.aov)
## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 platinum_type 1 74 1.591 0.211 0.015
## 2 sample_type 1 74 10.900 0.001 * 0.043
## 3 platinum_type:sample_type 1 74 1.302 0.258 0.005
bwtrim(CD8_Tcells ~ platinum_type*sample_type, id = pt_number, data = py$log_dens3_ptt)
## Call:
## bwtrim(formula = CD8_Tcells ~ platinum_type * sample_type, id = pt_number,
## data = py$log_dens3_ptt)
##
## value df1 df2 p.value
## platinum_type 3.8513 1 7.9619 0.0855
## sample_type 13.2831 1 6.5957 0.0091
## platinum_type:sample_type 3.1792 1 6.5957 0.1204
py$log_dens3_cr %>%
group_by(sample_type, chemo_regimen) %>%
shapiro_test(Macrophages)
## # A tibble: 4 × 5
## sample_type chemo_regimen variable statistic p
## <chr> <chr> <chr> <dbl> <dbl>
## 1 Cystectomy MVAC Macrophages 0.969 0.499
## 2 Cystectomy gemcitabine_platinum Macrophages 0.876 0.0000984
## 3 TUR MVAC Macrophages 0.881 0.00254
## 4 TUR gemcitabine_platinum Macrophages 0.957 0.0695
py$log_dens3_ptt %>%
group_by(sample_type, platinum_type) %>%
shapiro_test(Macrophages)
## # A tibble: 4 × 5
## sample_type platinum_type variable statistic p
## <chr> <chr> <chr> <dbl> <dbl>
## 1 Cystectomy Carboplatin Macrophages 0.923 0.458
## 2 Cystectomy Cisplatin Macrophages 0.902 0.0000603
## 3 TUR Carboplatin Macrophages 0.933 0.541
## 4 TUR Cisplatin Macrophages 0.926 0.000588
py$log_dens3_cr %>%
group_by(sample_type) %>%
levene_test(Macrophages ~ chemo_regimen)
## Warning: There were 2 warnings in `mutate()`.
## The first warning was:
## ℹ In argument: `data = map(.data$data, .f, ...)`.
## Caused by warning in `leveneTest.default()`:
## ! group coerced to factor.
## ℹ Run `dplyr::last_dplyr_warnings()` to see the 1 remaining warning.
## # A tibble: 2 × 5
## sample_type df1 df2 statistic p
## <chr> <int> <int> <dbl> <dbl>
## 1 Cystectomy 1 78 0.929 0.338
## 2 TUR 1 78 0.0228 0.880
box_m(py$log_dens3_cr[, "Macrophages", drop = FALSE], py$log_dens3_cr$chemo_regimen)
## # A tibble: 1 × 4
## statistic p.value parameter method
## <dbl> <dbl> <dbl> <chr>
## 1 0.812 0.367 1 Box's M-test for Homogeneity of Covariance Matric…
bwtrim(Macrophages ~ chemo_regimen*sample_type, id = pt_number, data = py$log_dens3_cr)
## Call:
## bwtrim(formula = Macrophages ~ chemo_regimen * sample_type, id = pt_number,
## data = py$log_dens3_cr)
##
## value df1 df2 p.value
## chemo_regimen 2.8689 1 59.9896 0.0955
## sample_type 13.0202 1 55.4999 0.0007
## chemo_regimen:sample_type 2.2379 1 55.4999 0.1403
res.aov <- anova_test(
data = py$log_dens3_cr, dv = Macrophages, wid = pt_number,
between = chemo_regimen, within = sample_type
)
get_anova_table(res.aov)
## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 chemo_regimen 1 78 2.467 0.120000 0.021
## 2 sample_type 1 78 14.688 0.000255 * 0.056
## 3 chemo_regimen:sample_type 1 78 5.646 0.020000 * 0.022
f1.ld.f1(y = py$log_dens3_cr$Macrophages, time = py$log_dens3_cr$sample_type, group = py$log_dens3_cr$chemo_regimen, subject = py$log_dens3_cr$pt_number, time.order = c('TUR', 'Cystectomy'), group.order = c('MVAC', 'gemcitabine_platinum'), description=F)$ANOVA.test
## F1 LD F1 Model
## -----------------------
## Check that the order of the time and group levels are correct.
## Time level: TUR Cystectomy
## Group level: MVAC gemcitabine_platinum
## If the order is not correct, specify the correct order in time.order or group.order.
## Statistic df p-value
## Group 3.184968 1 0.0743184625
## Time 14.764237 1 0.0001218242
## Group:Time 5.140979 1 0.0233674643
f1.ld.f1(y = py$log_dens3_ptt$Macrophages, time = py$log_dens3_ptt$sample_type, group = py$log_dens3_ptt$platinum_type, subject = py$log_dens3_ptt$pt_number, time.order = c('TUR', 'Cystectomy'), group.order = c('Cisplatin', 'Carboplatin'), description=F)$ANOVA.test
## F1 LD F1 Model
## -----------------------
## Check that the order of the time and group levels are correct.
## Time level: TUR Cystectomy
## Group level: Cisplatin Carboplatin
## If the order is not correct, specify the correct order in time.order or group.order.
## Statistic df p-value
## Group 9.6853793 1 0.001857400
## Time 7.2927779 1 0.006923234
## Group:Time 0.2585108 1 0.611145107
ps_ssgsea_cr = ps_ssgsea0[ps_ssgsea0.chemo_regimen.isin(['MVAC', 'gemcitabine_platinum'])]#.drop(['15 TUR', '15 Cystectomy', '7 TUR', '7 Cystectomy', '31 TUR', '31 Cystectomy'])
ps_ssgsea_ptt = ps_ssgsea0[ps_ssgsea0.platinum_type.isin(['Cisplatin', 'Carboplatin'])]
py$ps_ssgsea_cr %>%
group_by(sample_type, chemo_regimen) %>%
shapiro_test(TGFB_Mariathasan)
## # A tibble: 4 × 5
## sample_type chemo_regimen variable statistic p
## <chr> <chr> <chr> <dbl> <dbl>
## 1 Cystectomy MVAC TGFB_Mariathasan 0.936 0.164
## 2 Cystectomy gemcitabine_platinum TGFB_Mariathasan 0.964 0.144
## 3 TUR MVAC TGFB_Mariathasan 0.926 0.102
## 4 TUR gemcitabine_platinum TGFB_Mariathasan 0.989 0.921
py$ps_ssgsea_ptt %>%
group_by(sample_type, platinum_type) %>%
shapiro_test(TGFB_Mariathasan)
## # A tibble: 4 × 5
## sample_type platinum_type variable statistic p
## <chr> <chr> <chr> <dbl> <dbl>
## 1 Cystectomy Carboplatin TGFB_Mariathasan 0.889 0.134
## 2 Cystectomy Cisplatin TGFB_Mariathasan 0.957 0.0427
## 3 TUR Carboplatin TGFB_Mariathasan 0.896 0.167
## 4 TUR Cisplatin TGFB_Mariathasan 0.977 0.361
py$ps_ssgsea_cr %>%
group_by(sample_type) %>%
levene_test(TGFB_Mariathasan ~ chemo_regimen)
## Warning: There were 2 warnings in `mutate()`.
## The first warning was:
## ℹ In argument: `data = map(.data$data, .f, ...)`.
## Caused by warning in `leveneTest.default()`:
## ! group coerced to factor.
## ℹ Run `dplyr::last_dplyr_warnings()` to see the 1 remaining warning.
## # A tibble: 2 × 5
## sample_type df1 df2 statistic p
## <chr> <int> <int> <dbl> <dbl>
## 1 Cystectomy 1 68 0.454 0.503
## 2 TUR 1 68 6.11 0.0159
py$ps_ssgsea_ptt %>%
group_by(sample_type) %>%
levene_test(TGFB_Mariathasan ~ platinum_type)
## Warning: There were 2 warnings in `mutate()`.
## The first warning was:
## ℹ In argument: `data = map(.data$data, .f, ...)`.
## Caused by warning in `leveneTest.default()`:
## ! group coerced to factor.
## ℹ Run `dplyr::last_dplyr_warnings()` to see the 1 remaining warning.
## # A tibble: 2 × 5
## sample_type df1 df2 statistic p
## <chr> <int> <int> <dbl> <dbl>
## 1 Cystectomy 1 66 0.777 0.381
## 2 TUR 1 66 0.0422 0.838
box_m(py$ps_ssgsea_cr[, "TGFB_Mariathasan", drop = FALSE], py$ps_ssgsea_cr$chemo_regimen)
## # A tibble: 1 × 4
## statistic p.value parameter method
## <dbl> <dbl> <dbl> <chr>
## 1 3.47 0.0627 1 Box's M-test for Homogeneity of Covariance Matric…
box_m(py$ps_ssgsea_ptt[, "TGFB_Mariathasan", drop = FALSE], py$ps_ssgsea_ptt$chemo_regimen)
## # A tibble: 1 × 4
## statistic p.value parameter method
## <dbl> <dbl> <dbl> <chr>
## 1 3.36 0.0668 1 Box's M-test for Homogeneity of Covariance Matric…
res.aov <- anova_test(
data = py$ps_ssgsea_cr, dv = TGFB_Mariathasan, wid = pt_number,
between = chemo_regimen, within = sample_type
)
get_anova_table(res.aov)
## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 chemo_regimen 1 68 1.101 2.98e-01 0.010
## 2 sample_type 1 68 24.729 4.74e-06 * 0.123
## 3 chemo_regimen:sample_type 1 68 4.629 3.50e-02 * 0.026
tgf_df = bwtrim(TGFB_Mariathasan ~ chemo_regimen*sample_type, id = pt_number, data = py$ps_ssgsea_cr)
bwtrim(TGFB_Mariathasan ~ platinum_type*sample_type, id = pt_number, data = py$ps_ssgsea_ptt)
## Call:
## bwtrim(formula = TGFB_Mariathasan ~ platinum_type * sample_type,
## id = pt_number, data = py$ps_ssgsea_ptt)
##
## value df1 df2 p.value
## platinum_type 0.1571 1 8.7079 0.7014
## sample_type 26.5357 1 9.1856 0.0006
## platinum_type:sample_type 0.9076 1 9.1856 0.3652
f1.ld.f1(y = py$ps_ssgsea_cr$TGFB_Mariathasan, time = py$ps_ssgsea_cr$sample_type, group = py$ps_ssgsea_cr$chemo_regimen, subject = py$ps_ssgsea_cr$pt_number, time.order = c('TUR', 'Cystectomy'), group.order = c('MVAC', 'gemcitabine_platinum'), description=F)$ANOVA.test
## F1 LD F1 Model
## -----------------------
## Check that the order of the time and group levels are correct.
## Time level: TUR Cystectomy
## Group level: MVAC gemcitabine_platinum
## If the order is not correct, specify the correct order in time.order or group.order.
## Statistic df p-value
## Group 0.4801879 1 4.883372e-01
## Time 27.8105919 1 1.337916e-07
## Group:Time 3.1356103 1 7.659971e-02
f1.ld.f1(y = py$ps_ssgsea_ptt$TGFB_Mariathasan, time = py$ps_ssgsea_ptt$sample_type, group = py$ps_ssgsea_ptt$platinum_type, subject = py$ps_ssgsea_ptt$pt_number, time.order = c('TUR', 'Cystectomy'), group.order = c('Cisplatin', 'Carboplatin'), description=F)$ANOVA.test
## F1 LD F1 Model
## -----------------------
## Check that the order of the time and group levels are correct.
## Time level: TUR Cystectomy
## Group level: Cisplatin Carboplatin
## If the order is not correct, specify the correct order in time.order or group.order.
## Statistic df p-value
## Group 0.1589382 1 6.901360e-01
## Time 34.3513238 1 4.600873e-09
## Group:Time 1.7966710 1 1.801155e-01
pdl_cr = pdl1_df3_log[pdl1_df3_log.chemo_regimen.isin(['MVAC', 'gemcitabine_platinum'])]
pdl_ptt = pdl1_df3_log[pdl1_df3_log.platinum_type.isin(['Cisplatin', 'Carboplatin'])]
py$pdl_cr %>%
group_by(sample_type, chemo_regimen) %>%
shapiro_test(IC)
## # A tibble: 4 × 5
## sample_type chemo_regimen variable statistic p
## <chr> <chr> <chr> <dbl> <dbl>
## 1 Cystectomy MVAC IC 0.797 0.0000351
## 2 Cystectomy gemcitabine_platinum IC 0.883 0.000454
## 3 TUR MVAC IC 0.867 0.000994
## 4 TUR gemcitabine_platinum IC 0.904 0.00184
py$pdl_ptt %>%
group_by(sample_type, platinum_type) %>%
shapiro_test(IC)
## # A tibble: 4 × 5
## sample_type platinum_type variable statistic p
## <chr> <chr> <chr> <dbl> <dbl>
## 1 Cystectomy Carboplatin IC 0.932 0.570
## 2 Cystectomy Cisplatin IC 0.814 0.000000182
## 3 TUR Carboplatin IC 0.824 0.0698
## 4 TUR Cisplatin IC 0.895 0.0000616
py$pdl_cr %>%
group_by(sample_type) %>%
levene_test(IC ~ chemo_regimen)
## Warning: There were 2 warnings in `mutate()`.
## The first warning was:
## ℹ In argument: `data = map(.data$data, .f, ...)`.
## Caused by warning in `leveneTest.default()`:
## ! group coerced to factor.
## ℹ Run `dplyr::last_dplyr_warnings()` to see the 1 remaining warning.
## # A tibble: 2 × 5
## sample_type df1 df2 statistic p
## <chr> <int> <int> <dbl> <dbl>
## 1 Cystectomy 1 72 2.36 0.129
## 2 TUR 1 72 0.0119 0.913
box_m(py$pdl_cr[, "IC", drop = FALSE], py$pdl_cr$chemo_regimen)
## # A tibble: 1 × 4
## statistic p.value parameter method
## <dbl> <dbl> <dbl> <chr>
## 1 1.34 0.246 1 Box's M-test for Homogeneity of Covariance Matric…
#bwtrim(IC ~ chemo_regimen*sample_type, id = pt_number, data = py$pdl_cr)
#bwtrim(IC ~ chemo_regimen*sample_type, id = pt_number, data = py$pdl_ptt)
f1.ld.f1(y = py$pdl_cr$IC, time = py$pdl_cr$sample_type, group = py$pdl_cr$chemo_regimen, subject = py$pdl_cr$pt_number, time.order = c('TUR', 'Cystectomy'), group.order = c('MVAC', 'gemcitabine_platinum'), description=F)$ANOVA.test
## F1 LD F1 Model
## -----------------------
## Check that the order of the time and group levels are correct.
## Time level: TUR Cystectomy
## Group level: MVAC gemcitabine_platinum
## If the order is not correct, specify the correct order in time.order or group.order.
## Statistic df p-value
## Group 6.125007 1 1.332828e-02
## Time 19.517389 1 9.968807e-06
## Group:Time 2.295293 1 1.297667e-01
f1.ld.f1(y = py$pdl_ptt$IC, time = py$pdl_ptt$sample_type, group = py$pdl_ptt$platinum_type, subject = py$pdl_ptt$pt_number, time.order = c('TUR', 'Cystectomy'), group.order = c('Cisplatin', 'Carboplatin'), description=F)$ANOVA.test
## F1 LD F1 Model
## -----------------------
## Check that the order of the time and group levels are correct.
## Time level: TUR Cystectomy
## Group level: Cisplatin Carboplatin
## If the order is not correct, specify the correct order in time.order or group.order.
## Statistic df p-value
## Group 0.06227595 1 8.029342e-01
## Time 22.36575350 1 2.253569e-06
## Group:Time 0.72560845 1 3.943104e-01
res.aov <- anova_test(
data = py$pdl_cr, dv = IC, wid = pt_number,
between = chemo_regimen, within = sample_type
)
get_anova_table(res.aov)
## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 chemo_regimen 1 72 7.410 8.00e-03 * 0.056
## 2 sample_type 1 72 22.106 1.21e-05 * 0.116
## 3 chemo_regimen:sample_type 1 72 1.398 2.41e-01 0.008
cys_den = pdl1_df3_log[pdl1_df3_log.sample_type == 'Cystectomy'].set_index('pt_number')
tur_den = pdl1_df3_log[pdl1_df3_log.sample_type == 'TUR'].set_index('pt_number')
common_in = cys_den.index.intersection(tur_den.index)
delta_den = cys_den.loc[common_in].drop(['sample_type', 'chemo_regimen', 'platinum_type'], axis=1) - tur_den.loc[common_in].drop(['sample_type', 'chemo_regimen', 'platinum_type'], axis=1)
delta_den = delta_den.assign(pt_number = delta_den.index)
delta_den = delta_den.assign(chemo_regimen = delta_den.pt_number.map(pt_reg)).assign(platinum_type = delta_den.pt_number.map(platinum_type_dict))
log_dist3_cr = all_med_log[all_med_log.chemo_regimen.isin(['MVAC', 'gemcitabine_platinum'])]
log_dist3_ptt = all_med_log[all_med_log.platinum_type.isin(['Cisplatin', 'Carboplatin'])]
py$log_dist3_cr %>%
group_by(sample_type, chemo_regimen) %>%
shapiro_test(CD8_Tcell_to_tumor_median)
## # A tibble: 4 × 5
## sample_type chemo_regimen variable statistic p
## <chr> <chr> <chr> <dbl> <dbl>
## 1 Cystectomy MVAC CD8_Tcell_to_tumor_median 0.946 0.120
## 2 Cystectomy gemcitabine_platinum CD8_Tcell_to_tumor_median 0.963 0.122
## 3 TUR MVAC CD8_Tcell_to_tumor_median 0.779 0.0000217
## 4 TUR gemcitabine_platinum CD8_Tcell_to_tumor_median 0.951 0.0417
py$log_dist3_ptt %>%
group_by(sample_type, platinum_type) %>%
shapiro_test(CD8_Tcell_to_tumor_median)
## # A tibble: 4 × 5
## sample_type platinum_type variable statistic p
## <chr> <chr> <chr> <dbl> <dbl>
## 1 Cystectomy Carboplatin CD8_Tcell_to_tumor_median 0.936 0.574
## 2 Cystectomy Cisplatin CD8_Tcell_to_tumor_median 0.956 0.0168
## 3 TUR Carboplatin CD8_Tcell_to_tumor_median 0.912 0.366
## 4 TUR Cisplatin CD8_Tcell_to_tumor_median 0.839 0.000000407
py$log_dist3_cr %>%
group_by(sample_type) %>%
levene_test(CD8_Tcell_to_tumor_median ~ chemo_regimen)
## Warning: There were 2 warnings in `mutate()`.
## The first warning was:
## ℹ In argument: `data = map(.data$data, .f, ...)`.
## Caused by warning in `leveneTest.default()`:
## ! group coerced to factor.
## ℹ Run `dplyr::last_dplyr_warnings()` to see the 1 remaining warning.
## # A tibble: 2 × 5
## sample_type df1 df2 statistic p
## <chr> <int> <int> <dbl> <dbl>
## 1 Cystectomy 1 78 0.0740 0.786
## 2 TUR 1 78 1.78 0.187
box_m(py$log_dist3_cr[, "CD8_Tcell_to_tumor_median", drop = FALSE], py$log_dist3_cr$chemo_regimen)
## # A tibble: 1 × 4
## statistic p.value parameter method
## <dbl> <dbl> <dbl> <chr>
## 1 3.85 0.0497 1 Box's M-test for Homogeneity of Covariance Matric…
f1.ld.f1(y = py$log_dist3_cr$CD8_Tcell_to_tumor_median, time = py$log_dist3_cr$sample_type, group = py$log_dist3_cr$chemo_regimen, subject = py$log_dist3_cr$pt_number, time.order = c('TUR', 'Cystectomy'), group.order = c('MVAC', 'gemcitabine_platinum'), description=F)$ANOVA.test
## F1 LD F1 Model
## -----------------------
## Check that the order of the time and group levels are correct.
## Time level: TUR Cystectomy
## Group level: MVAC gemcitabine_platinum
## If the order is not correct, specify the correct order in time.order or group.order.
## Statistic df p-value
## Group 0.04884087 1 8.250927e-01
## Time 26.11852161 1 3.210892e-07
## Group:Time 0.04143483 1 8.387009e-01
f1.ld.f1(y = py$log_dist3_ptt$CD8_Tcell_to_tumor_median, time = py$log_dist3_ptt$sample_type, group = py$log_dist3_ptt$platinum_type, subject = py$log_dist3_ptt$pt_number, time.order = c('TUR', 'Cystectomy'), group.order = c('Cisplatin', 'Carboplatin'), description=F)$ANOVA.test
## F1 LD F1 Model
## -----------------------
## Check that the order of the time and group levels are correct.
## Time level: TUR Cystectomy
## Group level: Cisplatin Carboplatin
## If the order is not correct, specify the correct order in time.order or group.order.
## Statistic df p-value
## Group 0.04047519 1 8.405543e-01
## Time 22.55412414 1 2.043055e-06
## Group:Time 0.45706426 1 4.989992e-01
res.aov <- anova_test(
data = py$log_dist3_cr, dv = CD8_Tcell_to_tumor_median, wid = pt_number,
between = chemo_regimen, within = sample_type
)
get_anova_table(res.aov)
## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 chemo_regimen 1 78 0.035 8.52e-01 0.000323
## 2 sample_type 1 78 17.819 6.51e-05 * 0.061000
## 3 chemo_regimen:sample_type 1 78 0.212 6.47e-01 0.000769
py$log_dist3_cr %>%
group_by(sample_type, chemo_regimen) %>%
shapiro_test(Macrophage_to_tumor_median)
## # A tibble: 4 × 5
## sample_type chemo_regimen variable statistic p
## <chr> <chr> <chr> <dbl> <dbl>
## 1 Cystectomy MVAC Macrophage_to_tumor_median 0.846 4.10e-4
## 2 Cystectomy gemcitabine_platinum Macrophage_to_tumor_median 0.931 6.47e-3
## 3 TUR MVAC Macrophage_to_tumor_median 0.804 6.12e-5
## 4 TUR gemcitabine_platinum Macrophage_to_tumor_median 0.855 2.54e-5
py$log_dist3_ptt %>%
group_by(sample_type, platinum_type) %>%
shapiro_test(Macrophage_to_tumor_median)
## # A tibble: 4 × 5
## sample_type platinum_type variable statistic p
## <chr> <chr> <chr> <dbl> <dbl>
## 1 Cystectomy Carboplatin Macrophage_to_tumor_median 0.959 0.796
## 2 Cystectomy Cisplatin Macrophage_to_tumor_median 0.908 0.0000983
## 3 TUR Carboplatin Macrophage_to_tumor_median 0.789 0.0220
## 4 TUR Cisplatin Macrophage_to_tumor_median 0.804 0.0000000437
py$log_dist3_cr %>%
group_by(sample_type) %>%
levene_test(Macrophage_to_tumor_median ~ chemo_regimen)
## Warning: There were 2 warnings in `mutate()`.
## The first warning was:
## ℹ In argument: `data = map(.data$data, .f, ...)`.
## Caused by warning in `leveneTest.default()`:
## ! group coerced to factor.
## ℹ Run `dplyr::last_dplyr_warnings()` to see the 1 remaining warning.
## # A tibble: 2 × 5
## sample_type df1 df2 statistic p
## <chr> <int> <int> <dbl> <dbl>
## 1 Cystectomy 1 78 0.110 0.741
## 2 TUR 1 78 2.13 0.149
box_m(py$log_dist3_cr[, "Macrophage_to_tumor_median", drop = FALSE], py$log_dist3_cr$chemo_regimen)
## # A tibble: 1 × 4
## statistic p.value parameter method
## <dbl> <dbl> <dbl> <chr>
## 1 4.33 0.0375 1 Box's M-test for Homogeneity of Covariance Matric…
f1.ld.f1(y = py$log_dist3_cr$Macrophage_to_tumor_median, time = py$log_dist3_cr$sample_type, group = py$log_dist3_cr$chemo_regimen, subject = py$log_dist3_cr$pt_number, time.order = c('TUR', 'Cystectomy'), group.order = c('MVAC', 'gemcitabine_platinum'), description=F)$ANOVA.test
## F1 LD F1 Model
## -----------------------
## Check that the order of the time and group levels are correct.
## Time level: TUR Cystectomy
## Group level: MVAC gemcitabine_platinum
## If the order is not correct, specify the correct order in time.order or group.order.
## Statistic df p-value
## Group 0.8742298 1 3.497870e-01
## Time 23.2285838 1 1.438434e-06
## Group:Time 0.1916301 1 6.615633e-01
f1.ld.f1(y = py$log_dist3_ptt$Macrophage_to_tumor_median, time = py$log_dist3_ptt$sample_type, group = py$log_dist3_ptt$platinum_type, subject = py$log_dist3_ptt$pt_number, time.order = c('TUR', 'Cystectomy'), group.order = c('Cisplatin', 'Carboplatin'), description=F)$ANOVA.test
## F1 LD F1 Model
## -----------------------
## Check that the order of the time and group levels are correct.
## Time level: TUR Cystectomy
## Group level: Cisplatin Carboplatin
## If the order is not correct, specify the correct order in time.order or group.order.
## Statistic df p-value
## Group 0.07567483 1 7.832466e-01
## Time 35.48936404 1 2.564470e-09
## Group:Time 2.54512061 1 1.106355e-01
res.aov <- anova_test(
data = py$log_dist3_cr, dv = Macrophage_to_tumor_median, wid = pt_number,
between = chemo_regimen, within = sample_type
)
get_anova_table(res.aov)
## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 chemo_regimen 1 78 1.185 2.80e-01 0.010000
## 2 sample_type 1 78 18.644 4.58e-05 * 0.080000
## 3 chemo_regimen:sample_type 1 78 0.177 6.75e-01 0.000828